1. INTRODUCTION

Este trabajo es una continuación de este, por lo tanto lo primero que tendremos que tener claro, para poder interpretar los diferentes clusters que vayamos sacando y tener un criterio aparte de los matemáticos para quedarnos con uno u otro tipo de cluster, será recordar cuál era el objetivo para el que se recogieron las variables que se utilizaron en esa primera parte del proyecto y que se utilizarán también en esta. El objetivo era, en definitiva, ver cómo se posicionaban los países en función del spread y el bounty.

Históricamente se ha hablado siempre del primer, segundo y tercer mundo; sin embargo, en los últimos años se habla más de países desarrollados (bien posicionados respecto al efecto de spread generado por el progreso tecnológico) y países en vías de desarrollo (los que se están quedando atrás con respecto al resto debido a este progreso tecnológico, ya que se están comiendo una parte más pequeña de la tarta que este genera).

Por lo tanto, estaremos cómodos tanto con 2 como con 3 clusters. Es posible que encontremos alguna composición de un número mayor de clusters con la que nos sintamos cómodos porque al analizarla tenga bastante sentido y una interpretabilidad suficiente. Sin embargo, a priori buscaremos encontrar 2 o 3 clusters.

Como todos los gráficos de cluster los vamos a realizar sobre los dos primeros componentes principales, remito al lector a este link, en el que se especifica y se analizan los primeros componentes principales.

El mejor cluster de cada tipo (Kmeans, Hierarchical, etc) se evaluará también sobre otros subconjuntos de variables dentro de la muestra de 21 variables con la que vamos a trabajar.

Adicionalmente, incluimos un pequeño párrafo con las variables que son objeto de estudio:

  1. Gdp Per Capita: la sacamos de WB
  2. Inflation: la sacamos de WB. Medida como índice de precios al consumo.
  3. Gni Per Capita: la sacamos de WB
  4. Urban Population: porcentaje de población urbana. La sacamos de UN.
  5. Employment Agriculture: porcentaje de empleo en agricultura. UN
  6. Economy Agriculture: Importancia de la agricultura en la economía. UN.
  7. Urban Population Growth: de 2016 a 2017. UN.
  8. Credit Information: Del 1 al 8 la información de crédito que tiene la gente de un país. WB.
  9. Productivity: GDP/Hour worked. WB.
  10. Renewable Energy: Uso de energía renovable. WB.
  11. Renewable Electricity: Porcentaje de la energía eléctrica que es renovable. WB.
  12. Electricity 2016: WB.
  13. Clean Fuels 2016: Acceso a “clean fuels”. WB.
  14. Infant Mortality Rate: WHO.
  15. Sanitation Services: WHO.
  16. Deaths for Household Pollution: WHO.
  17. Water: WHO.
  18. Fertility Rate: UN.
  19. Account Ownership: WB.
  20. Life Expectancy: UN.
  21. Vulnerable Employment: WB.

Eliminamos la primera columna, ya que sólo contiene el índice. También eliminamos la observación en la que la inflación era máxima, pues, como se vio en el tratamiento previo de las variables, que se puede encontrar aqui, esta observación era un outlier.

En este caso no realizamos un análisis exploratorio de las variables, porque eso ya lo hicimos en la preparación (ver link que aparece arriba); por lo tanto iremos directamente a probar diferentes métodos de cluster, y buscaremos el que más nos convenza a la hora de separar los países en estos diferentes grupos.

Nos quedamos sólo con las columnas numéricas para la matriz X.

2. Visualización Multivariante.

Antes de entrar a implementar las técnicas de aprendizaje no supervisado con las que trataremos de agrupar los países, vamos a visualizar de forma sencilla la matriz de correlaciones de las variables, así como los histogramas de la distribución de estas variables. Esto último es lo primero que vamos a hacer, pues así podremos ver qué transformaciones en las variables podríamos hacer para tratar de asemejarlas más a una normal. En este link vemos formas típicas de conseguir que distribuciones asimétricas tomen una forma más parecida a la normal. Además de las que vienen en el enlace, probaremos otras para quedarnos con las que más nos guste, siempre usando el histograma de la variable transformada como criterio. Sine embargo, con el fin de no ensuciar el report, mostraremos directamente sólo el código de transformación de estas variables, sin dibujar de nuevo estos histogramas.

En este chunk se pueden ver las transformaciones que les hacemos a las variables para tratar de que sus distribuciones se parezcan lo máximo posible a la Normal. Esto será importante especialmente a la hora de interpretar las probabilidades de los Mixture Models, que probaremos más adelante.

Transformaciones.

X[ , 1] <- X[ , 1]**(1/3)

#X[ , 2] <- X[ , 2]**(1/3)

X[ , 3] <- X[ , 3]**(1/4)

X[ , 5] <- X[ , 5]**(1/3)

X[ , 8] <- X[ , 8]**(1/2)

X[,9] <- X[ , 9]**(2)**(1/3)

X[, 11] <- X[,11]**(1/3)

X[ , 13] <- X[,13]**(1/2)

X[, 15] <- X[,15]**(1/4)

X[,16] <- X[,16]**(1/2)

#X[,17] <- X[,17]**(1/4)

X[,19] <- X[,19]**(1/3)


X <- data.frame(scale(X))

Podemos ver que la mayoría de las variables tienen bastante correlación, tanto positiva como negativa, con el resto. Esto es buena señal, pues indica que, pares a pares (es decir, en dimensión 2), nuestras variables están bastante relacionadas; esto hará que los modelos de aprendizaje no supervisado que implementaremos a continuación tengan mayor facilidad para separar a los países en grupos.

3. K-MEANS CLUSTERING

Using Euclidean Distance

Probamos primero con 5 clusters, pese a que, como ya se comentó en la introducción, buscamos preferiblemente 3 clusters (o entre 2 y 3).

Para cada método de cluster, se enseñará sólo el código de generación para el primer tipo que se pruebe de ese método; en los demás es igual pero cambiando los parámetros como se indicará en los títulos de los apartados o en el propio texto descriptivo de los resultados o nuevas pruebas.

fit = kmeans(X, centers=5, nstart=100)
groups = fit$cluster
#groups
barplot(table(groups), col="blue")

centers=fit$centers
centers
##   gdp_per_cap  inflation gni_per_cap    urban_p productivity  credit_inf
## 1 -0.32684608 -0.2665991 -0.28436011 -0.1384272   -0.2125053  0.47140805
## 2  1.40617124 -0.5015226  1.38975542  1.0172553    1.2224013  0.54036621
## 3  0.02419864  0.3327386  0.08761484  0.2245190    0.3409200  0.01201352
## 4 -1.11412294  0.9075402 -1.18388174 -0.9423741   -1.3117730 -0.94862434
## 5 -0.74670445 -0.1864966 -0.79203840 -0.8211653   -0.8148517 -0.61891526
##      renew_en renew_electr electricity_2016 clean_fuels_2016 inf_mort_rate
## 1 -0.02238447   0.44398477        0.5327119        0.1667835    -0.1420763
## 2 -0.48843966  -0.04074668        0.6239364        0.8522864    -1.1611765
## 3 -0.87348254  -0.78239308        0.5818854        0.7360304    -0.1981290
## 4  1.59717242   0.32985108       -1.8265185       -1.5472748     1.3701156
## 5  0.30659661   0.10239740       -0.6629780       -1.0263136     0.9862122
##   sanitation_services deaths_household_pollution      water  fert_rate
## 1           0.3162657                  0.2922737  0.3766224 -0.2978457
## 2           0.8436136                 -1.1339905  0.7753450 -0.9189310
## 3           0.5867365                 -0.5722864  0.4844916 -0.2879624
## 4          -1.6589471                  1.1896931 -1.7622550  1.5181949
## 5          -0.9095397                  0.9925236 -0.6274547  0.7596367
##     life_exp urban_pop_growth      emp_agr    eco_agr vulner_empl
## 1  0.2542822       -0.2868240  0.002811144  0.2569874  0.04609847
## 2  1.0281593       -0.6359001 -0.948799884 -1.0873899 -1.01792821
## 3  0.2473097       -0.3302692 -0.539902591 -0.4052800 -0.52012687
## 4 -1.5232020        1.2979051  1.517756967  1.3075144  1.54357537
## 5 -0.8223502        0.5529592  0.700815723  0.6052136  0.69963044
##       acc_own
## 1 -0.28961642
## 2  1.21761600
## 3  0.06864725
## 4 -0.91305480
## 5 -0.75933271
# clusplot
fviz_cluster(fit, data = X, geom = c("point"),ellipse.type = 'norm', pointsize=1)+
  theme_minimal()+geom_text(label=countries,hjust=0, vjust=0,size=2,check_overlap = T)+scale_fill_brewer(palette="Paired")

# Silhouette plot
# The silhouette value in [-1,1] measures the similarity (cohesion) of a data point to its cluster relative to other clusters (separation). 
# Silhouette plots rely on a distance metric and suggest that the data matches its own cluster well.
# The larger the silhouette widths, the better.
d <- dist(X, method="euclidean")  
sil = silhouette(groups, d)
plot(sil, col=1:5, main="", border=NA)

summary(sil)
## Silhouette of 182 units in 5 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##         41         47         37         30         27 
## 0.20516587 0.39670270 0.04886944 0.19057532 0.13933644 
## Individual silhouette widths:
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -0.08574  0.09069  0.19322  0.21068  0.30593  0.54489

Vemos que el average silhouette width no es demasiado alto; probaremos bajando el número de clusters.

Además, nos damos cuenta en este punto de que las transformaciones que hemos realizado previamente de las variables han “estropeado” los datos, en el sentido de que no queda tan clara la división, basada en los dos primeros componentes principales, como cuando utilizábamos los datos sin transformar. Por lo tanto, recogemos cuerda y utilizaremos los datos sin transformar, pues con estos los resultados son bastante más interpretables. Eso puede suceder porque, pese a que las variables por separado se comporten de forma más parecida a la normal, la distribución multidimensional de todo el dataset no lo haga.

Bajamos ahora el número de clusters a 4.

## Silhouette of 182 units in 4 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        36        29        33        84 
## 0.3919825 0.1839574 0.1496398 0.2296468 
## Individual silhouette widths:
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -0.07903  0.14221  0.23234  0.23997  0.32905  0.54965

Sube un poco el average silhouette width. Seguimos bajando en el número de centroides para ver si mejora. Probaremos ahora con 3 centroides.

## Silhouette of 182 units in 3 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        39        55        88 
## 0.3987943 0.3057256 0.2464965 
## Individual silhouette widths:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -0.0734  0.1934  0.3212  0.2970  0.4099  0.5617

Vemos que con 3 clusters el average silhouette width sigue mejorando; además, estos clusters, visualmente, por los países que forman los clusters, podríamos decir que representarían el primer, segundo y tercer mundo respectivamente, aunque antes de saltar a las conclusiones debemos hacer un estudio detallado de las características de cada cluster para poder afirmar esto.

##   gdp_per_cap  inflation gni_per_cap     urban_p productivity credit_inf
## 1   1.6211935 -0.5446340   1.6543070  1.10617407    1.4874171  0.5075753
## 2  -0.6356268  0.4246104  -0.6539681 -0.87108893   -0.8400000 -0.8068516
## 3  -0.3212167 -0.0240096  -0.3244287  0.05419435   -0.1341962  0.2793341
##     renew_en renew_electr electricity_2016 clean_fuels_2016 inf_mort_rate
## 1 -0.5579749   -0.1685010        0.5993024        0.8531021    -0.9127375
## 2  1.0390537    0.3440539       -1.2998611       -1.3107411     1.2595349
## 3 -0.4021242   -0.1403571        0.5468133        0.4411339    -0.3827006
##   sanitation_services deaths_household_pollution      water  fert_rate
## 1           0.8557941                 -0.9417420  0.7846983 -0.8116898
## 2          -1.3447276                  1.1469979 -1.2734136  1.2176491
## 3           0.4611823                 -0.2995108  0.4481194 -0.4013045
##     life_exp urban_pop_growth    emp_agr    eco_agr vulner_empl
## 1  1.1401976       -0.5241457 -0.9815974 -0.8516165  -1.0470265
## 2 -1.2200850        0.9950285  1.1614203  1.0139171   1.1546237
## 3  0.2572383       -0.3896010 -0.2908616 -0.2562773  -0.2576167
##       acc_own
## 1  1.30387480
## 2 -0.85861219
## 3 -0.04122098

Como vemos, el primer cluster tiene un gdp per capita bastante inferior ala media, una inflación bastante superior, gni per capita muy bajo, población urbana, productividad e información al crédito aún más bajas. Por otro lado, el uso de energías renovables es muy alto (siempre con respecto a la media), el acceso a electricidad renovable también alto aunque no tanto, el uso de electricidad y el acceso a clean fuels en 2016 ambos muy bajos, una alta tasa de mortalidad infantil (altísima), un acceso a los servicios sanitarios bajísimo con respecto a los demás, muertes por household pollution muy altas, bajísimo acceso a agua limpia, una alta tasa de fertilidad, muy baja esperanza de vida, alto crecimiento de la población urbana, alto empleo y economía agrícola, un muy alto número de empleos vulnerables, y muy baja account ownership. Como podemos ver, este primer cluster correspondería al tercer mundo claramente. Serían aquellos países que peor posicionados están con respecto al spread.

Por otro lado, el segundo cluster sigue teniendo un gdp per capita bajo, una inflación cercana a la media, un gni per capita bajo (aunque no tan bajo como el de los países del cluster 1), una poblacion urbana normal (media), productividad cercana a la media, información al crédito alta, aunque no tan alta como en el cluster 3, bajo uso de energías renovables, bajo uso de electricidad renovable aunque cercano a la media. Por otro lado, estos países tienen alto acceso a la electricidad y a comnustibles limpios en 2016, una mortalidad infantil moderadamente baja, un acceso a servicios sanitarios moderadamente alto, bajas muertes por household pollution, moderadamente alto acceso al agua, una tasa de fertilidad baja, aunque no tan baja como la de los países del cluster 3, un crecimiento de la población urbana moderadamente bajo; tanto en empleo agrícola, como en agricultura agrícola, como en empleos vulnerables tienen una puntuación ligeramente inferior a la media, y un porcentaje de personas con cuenta bancaria muy cercano a la media. Este cluster correspondería a los países del segundo mundo, aquellos que están posicionados más o menos en la media con respecto al spread; ni se están aprovechando tanto de las consecuencias del progreso tecnológico como los del tercer clusters ni están padeciendo tanto estos efectos como los del primer cluster.

Por último, el tercer cluster tiene un gdp per capita, gni per capita, población urbana, productividad muy altas; una inflación bastante inferior a la media. Tiene la mayor información crediticia de todos los clusters, un bajo uso de energías renovables, un uso de electricidad renovable muy similar a la del segundo cluster, un acceso a la electricidad alto y acceso a los combustibles limpios similar al del cluster 2. La mortalidad infantil en este cluster es la más baja de todas, mientras que el acceso a servicios sanitarios es bastante alto, las muertes por household pollution es muy baja, el acceso al agua limpia es el más alto de todos los clusters, la tasa de fertilidad es muy baja, mientras que la esperanza de vida es altísima, muy por encima de los otros dos clusters. Tienen el crecimiento de la población urbana más bajo de todos los clusters, y el peso tanto del trabajo como de la economía agrícola es muy bajo en sus sistemas de producción. El porcentaje de empleo vulnerable es muchisimo más bajo que en el resto de clusters, y llama la atención el porcentaje de personas con cuenta bancaria, en cuyo apartado están claramente posicionados muy positivamente, con mucha diferencia con respecto a los otros dos clusters. Este cluster podría ser identificado, por lo tanto, como aquellos países que mejor posicionados están con respecto al efecto del spread; es decr, el primer mundo.

Estamos cómodos con esta composición de clusters, especialmente después de observar cuál es la media de los países de cada cluster en cada una de las variables introducidas, pues nos cuadra con la idea que teníamos al principio de que probablemente encontráramos tres grupos en este aspecto: aquellos que están quedándose claramente atrás debido al progreso tecnológico, los que están más o menos en la media, es decir que no sufren tanto las consecuencias de este progreso pero tampoco se están distanciando positivamente de los demás debido a este; y por último los países que están sacando una ventaja significativa del resto en la calidad de vida, es decir aquellos países que están muy bien posicionados con respecto al spread.

Vamos a probar con 2 centroides:

## Silhouette of 182 units in 2 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        59       123 
## 0.3541759 0.4486794 
## Individual silhouette widths:
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -0.01397  0.34020  0.45628  0.41804  0.52190  0.61582

Vemos que con 2 clusters parece que incluso mejora con respecto a 3.

##   gdp_per_cap  inflation gni_per_cap    urban_p productivity credit_inf
## 1  -0.6331667  0.3760846  -0.6501957 -0.8903473   -0.8330055 -0.7459775
## 2   0.3037141 -0.1803983   0.3118825  0.4270771    0.3995717  0.3578266
##     renew_en renew_electr electricity_2016 clean_fuels_2016 inf_mort_rate
## 1  0.9920984    0.3636492       -1.1917260       -1.2550132     1.1963608
## 2 -0.4758846   -0.1744333        0.5716409        0.6019982    -0.5738641
##   sanitation_services deaths_household_pollution     water  fert_rate
## 1          -1.2538326                  1.1193596 -1.183160  1.1334424
## 2           0.6014319                 -0.5369286  0.567532 -0.5436837
##     life_exp urban_pop_growth    emp_agr    eco_agr vulner_empl    acc_own
## 1 -1.1542185        0.9343703  1.1373573  0.9968732   1.1343759 -0.8069055
## 2  0.5536495       -0.4481939 -0.5455616 -0.4781750  -0.5441315  0.3870523

Sin embargo, pese a que matemáticamente parecen encajar mejor dos clusters que tres, al observar la media que tiene cada uno de los clusters en las diferentes variables nos encontraríamos más cómodos teniendo 3 grupos, pues por ejemplo en las variables económicas, en este caso, no se observa tanta diferencia (porque agrupa a los países del segundo mundo en el primero o en el tercero sin pertenecer del todo a ninguno de estos dos, lo cual hace que las diferencias en la media de estas variables no sea tan significativa). Mientras que con 3 clusters sí podíamos identificar una gran diferencia del tercer cluster (primer mundo) con respecto a los otros dos en algunas variables como el porcentaje de adultos con cuenta bancaria, o las variables económicas ya mencionadas antes, además de otras variables que miden la calidad de vida como el acceso al agua o a la sanidad, en este caso las diferencias son más sutiles porque los dos clusters están formados por países más dispares entre sí (Brasil y Luxemburgo, por ejemplo, están en el mismo cluster, y mientras que el primero está aún a años luz (aunque tenga mucho potencial) de posicionarse bien con respecto al spread, el segundo podríamos decir que es el país más desarrollado del mundo, o uno de ellos, junto a Singapur, Islandia…).

En caso de utilizar un k-means cluster, por lo tanto, ¿con cuántos nos quedamos?:

Ambos gráficos nos vienen a decir lo que acabábamos de ver nosotros probando y viendo el average silhouette width, es decir, que lo óptimo parecen ser 2 clusters. Esto, de todas formas, es sólo de manera matemática, y como dijimos antes nos sentiríamos más cómodos cogiendo 3 clusters.

Pruebas con Submuestras

Recordando lo que hicimos con los datos en el estudio de las variables y el ajuste de PCA y Factor Analysis, también podríamos eliminar algunas variables de nuestro set de datos que no tienen tanta relación con las demás y podrían estar incluyendo ruido en nuestra muestra. Esto es lo que se comentaba en la introducción, probaremos cada tipo de cluster con muestras de datos reducidas. Esta es la primera.

Con estos datos reducidos, cogiendo 3 centroides, quedaría asi:

## Silhouette of 182 units in 3 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        53        91        38 
## 0.3851991 0.3602663 0.4495979 
## Individual silhouette widths:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -0.0574  0.2872  0.4304  0.3862  0.5154  0.6165

Sin embargo, pese a observar que los países que forman cada cluster tendrían más en común entre sí, pues el average silhouette width sube, decidimos quedarnos con los datos originales, es decir con la matriz X. Es evidente que al reducir la dimensionalidad de esta matriz (quitando columnas), tendríamos mejores resultados, pero esto puede deberse únicamente al bajo número de observaciones que tenemos, y no tanto a la poca relevancia de esas variables. Por este motivo decidimos trabajar con la matriz completa aunque esto implique incluir un poco de ruido en nuestros datos, lo cual, como comentamos, es inevitable debido a que \(\frac{p}{n}\) es grande. Seguiremos probando con submuestras de datos, de todas formas, pues esto nos da información sobre cuáles son las variables que ayudan a agrupar mejor los datos de los países.

Cogemos ahora otra submuestra, sólo con las variables de carácter más económico. Además, dejamos preparada la tercera submuestra, que contendrá aquellas variables más relacionadas con la salud y la calidad de vida como tal, sin tener en cuenta lo económico.

Veamos el resultado del cluster k-means con 3 centroides para las variables económicas:

X_econ <- X[, c("gdp_per_cap", "inflation", "gni_per_cap", "productivity", "credit_inf", "emp_agr", "eco_agr", "vulner_empl", "acc_own")]

X_health <- X[, c("renew_en", "renew_electr", "electricity_2016", "clean_fuels_2016", "inf_mort_rate", "sanitation_services", "deaths_household_pollution", "water", "fert_rate", "life_exp")]

fit = kmeans(X_econ, centers=3, nstart=100)
groups = fit$cluster
#groups
barplot(table(groups), col="blue")

centers=fit$centers
#centers

# clusplot
fviz_cluster(fit, data = X_econ, geom = c("point"),ellipse.type = 'norm', pointsize=1)+
  theme_minimal()+geom_text(label=countries,hjust=0, vjust=0,size=2,check_overlap = T)+scale_fill_brewer(palette="Paired")

d <- dist(X_econ, method="euclidean")  
sil = silhouette(groups, d)
plot(sil, col=1:4, main="", border=NA)

summary(sil)
## Silhouette of 182 units in 3 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        89        60        33 
## 0.3415002 0.2513332 0.4933026 
## Individual silhouette widths:
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -0.06966  0.21879  0.36382  0.33930  0.45695  0.64764
i=1  # plottinng the centers in cluster 1
barplot(centers[i,], las=2, col="darkblue", main=paste("Cluster", i))

i=2  # plottinng the centers in cluster 1
barplot(centers[i,], las=2, col="darkblue", main=paste("Cluster", i))

i=3

barplot(centers[i,], las=2, col="darkblue", main=paste("Cluster", i))

require(rworldmap)

matched <- joinCountryData2Map(data.frame(country=countries, value=groups), joinCode="NAME", nameJoinColumn="country")
## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data
mapCountryData(matched, nameColumnToPlot="value", mapTitle="Clusters for Countries in the World by K-Means, only with economic variables", 
               catMethod = "categorical", colourPalette = "topo")

Vemos que los resultados no cambian en exceso con respecto a cuando usamos todas las variables. Esto nos indica la gran importancia de las variables económicas a la hora de distinguir a unos países de otros. Hemos decidido hacer 3 clusters porque nos parece más adecuado e interpretable dividir a los países en primer, segundo y tercer mundo que solo en primer y segundo mundo. Nos llama la atención, por ejemplo, que Portugar quedaría en el segundo mundo en términos meramente económicos (en España parece que no estamos tan mal y seguiríamos formando parte de los países más desarrollados también sólo en términos económicos). India caería al tercer mundo también según el clustering de K-means generado con estas variables solo económicas.

Ahora haremos lo mismo con las variables relacionadas con la salud y la calidad de vida en general (mortalidad infantil, tasa de fertilidad, esperanza de vida…):

## Silhouette of 182 units in 3 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        37       111        34 
## 0.1323641 0.5194394 0.3274172 
## Individual silhouette widths:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -0.1234  0.2541  0.4248  0.4049  0.6015  0.6667

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Curiosamente aquí si encontramos diferencias muy significativas con respecto a incluir todas las variables, pues hay muy pocos países (principalmente africanos y asiáticos) que pertenezcan al tercer mundo en términos de variables relacionadas con la salud y con la calidad de vida, mientras que la gran mayoría de países entran dentro del primer mundo. Esto viene a decir que las diferencias actuales entre los países en términos de salud y calidad de vida son bastante pequeñas, y nos muestra que las variables económicas son las de mayor peso a la hora de determinar el posicionamiento de los países con respecto al spread. Esto viene demostrado en el silhouette width del tercer cluster (primer mundo), ya que es de 0.52, bastante mayor a la de los otros dos clusters. 111 países de los 182 estudiados caen dentro de este.

Nota: el número de cluster puede variar con respecto a los comentarios, ya que cada vez que se ejecuta y se compila el algoritmo les asigna un número aleatorio a cada cluster. En el momento actual, el primer mundo sería el tercer cluster, el segundo mundo el primer cluster y el tercer mundo el segundo cluster; pero esto puede variar cada vez que se ejecute. Como tenemos claros algunos países que deben estar en cada uno de los grupos, guíese el lector por estos a la hora de interpretar los resultados, pues es posible que se aprecien contradicciones entre la redacción y los gráficos. Por dar algunas pistas: USA siempre forma parte del primer mundo en todos los modelos utilizados, los países del centro y sur de África forman parte generalmente del tercer mundo, y países como Mongolia forman parte generalmente del segundo. Remito al lector a estas referencias debido a que en caso contrario tendría que editar cada vez que ejecuto el texto todo lo escrito, lo cual, aparte de costoso temporalmente, seguramente provocaría errores e incongruencias.

Using Mahalanobis Distance

Hasta este punto hemos estado haciendo cluster con las distancias euclídeas, pero vamos a probar a hacerlo con la distancia de Mahalanobis. Para ello, realizamos las transformaciones a la matriz X como vemos en el chunk de código inferior. Esta transformación “truca” los datos de tal forma que al calcularse las distancias euclídeas de esta nueva matriz, es como si se estuviera calculando la distancia de Mahalanobis sobre los datos originales.

# k-means with Mahalanobis distance

S_x <- cov(X_unscaled)
iS <- solve(S_x) #para coger la inversa de la matrix
e <- eigen(iS)
V <- e$vectors
B <- V %*% diag(sqrt(e$values)) %*% t(V)
Xtil <- scale(X_unscaled,scale = FALSE)
XS <- Xtil %*% B

# kmeans with 2 clusters
fit = kmeans(XS, centers=2, nstart=100)
groups = fit$cluster
barplot(table(groups), col="blue")

centers=fit$centers
colnames(centers)=colnames(X)
#centers

# clusplot
fviz_cluster(fit, data = X, geom = c("point"),ellipse.type = 'norm', pointsize=1)+
  theme_minimal()+geom_text(label=countries,hjust=0, vjust=0,size=2,check_overlap = T)+scale_fill_brewer(palette="Paired")

En los gráficos de arriba queda suficientemente claro que realizar un k-means clusters utilizando la distancia de Mahalanobis nos da peores resultados que utilizando la distancia euclídea. Por este motivo no entraremos a explicar las características de cada uno de estos clusters, e iremos directamente a implementar el siguiente método de aprendizaje no supervisado: PAM; no sin antes visualizar los resultados del modelo de cluster por el método de k-means que más nos cuadra de los probados arriba.

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Como vemos en el mapa del mundo, nos cuadran bastante los países que pertenecen al primer mundo (cluster 3), segundo mundo (cluster 2) y tercer mundo (cluster 1). Probaremos ahora con el siguiente método: PAM.

4. PAM

Using Euclidean Distance

# number of clusters?
fviz_nbclust(X, pam, method = 'wss')

fviz_nbclust(X, pam, method = 'silhouette')

# PAM clustering, based on medioids
fit.pam <- pam(X, k=2, metric = "euclidean") # 2 cluster solution
fit.pam$medoids
##    gdp_per_cap  inflation gni_per_cap    urban_p productivity credit_inf
## 33  -0.6435012 -0.7046321  -0.6525326 -0.1063936   -0.8740631 -1.2085410
## 42   0.3844251 -0.3483195   0.3373432  0.6971845    0.6696516  0.7408338
##      renew_en renew_electr electricity_2016 clean_fuels_2016 inf_mort_rate
## 33  1.5459515    1.3198344       -0.7815197       -1.1109502     1.6683661
## 42 -0.6093636   -0.6762417        0.6553683        0.8557089    -0.9702128
##    sanitation_services deaths_household_pollution      water  fert_rate
## 33          -1.1885475                  1.3631735 -1.1869961  1.4674467
## 42           0.8566682                 -0.8092583  0.7934034 -0.9611414
##      life_exp urban_pop_growth   emp_agr    eco_agr vulner_empl    acc_own
## 33 -1.8001742         0.895438  1.535388  0.8394029   1.3416915 -0.9453004
## 42  0.8915702        -1.048183 -1.014676 -0.7953806  -0.8615597  0.8047001
#fit.pam$clustering


fviz_cluster(fit.pam, data = X, geom = c("point"), ellipse.type = 'norm', pointsize=1)+
  theme_minimal()+geom_text(label=countries,hjust=0, vjust=0,size=2,check_overlap = T)+scale_fill_brewer(palette="Paired")

Vemos que de nuevo, matemáticamente, parece que lo que más sentido tiene es tomar 2 clusters. Probaremos con otras distancias, así como con 3 clusters, y explicaremos los resultados de aquella composición que más nos encaje para el objetivo de este estudio.

##    gdp_per_cap  inflation gni_per_cap    urban_p productivity credit_inf
## 33  -0.6435012 -0.7046321  -0.6525326 -0.1063936   -0.8740631 -1.2085410
## 38  -0.5477250 -0.6763866  -0.5564911  0.3731611   -0.6336868  0.4159380
## 51   0.8075529 -0.4456493   0.8688122  0.9823251    1.1591846  0.7408338
##      renew_en renew_electr electricity_2016 clean_fuels_2016 inf_mort_rate
## 33  1.5459515   1.31983443       -0.7815197       -1.1109502     1.6683661
## 38 -0.1989357  -0.40468249        0.3894632        0.1634364    -0.3470056
## 51 -0.5595744   0.05001927        0.6553683        0.9318907    -0.9702128
##    sanitation_services deaths_household_pollution        water  fert_rate
## 33          -1.1885475                 1.36317353 -1.186996132  1.4674467
## 38          -0.3022874                -0.01661423  0.001243579 -0.2672591
## 51           0.8907551                -0.98540139  0.793403387 -1.0999179
##      life_exp urban_pop_growth     emp_agr    eco_agr vulner_empl
## 33 -1.8001742       0.89543799  1.53538810  0.8394029  1.34169151
## 38  0.1366442      -0.04692393  0.04749208 -0.2288715  0.02836333
## 51  1.4237312      -0.81259299 -0.95447259 -0.7872876 -0.94770847
##       acc_own
## 33 -0.9453004
## 38 -0.5864026
## 51  1.2860648

La composición con 3 medioides nos encaja mejor que con 2; probemos con otras distancias antes de decidir con cuál nos quedamos.

Using Manhattan Distance

Vemos que en el caso de PAM utilizando la distancia de Manhattan, lo óptimo (matemáticamente) también sería coger dos clusters. Con la distancia de Manhattan el average silhouette width sale algo más alto que con la distancia euclídea, sin importar el número de clusters (para todos en general es más alta).

Veamos los resultados para PAM, usando distancia de Manhattan y con 3 medioides.

Con 3 clusters nos ocurre lo mismo que comentábamos antes, y es que nos encajan más los resultados que cuando agrupamos a lo países solamente en 2 grupos. Vemos que utilizando el método de PAM con la distancia Manhattan mejoramos los resultados para este cluster de 3; midiéndolo como el average silhouette width y comparándolo con PAM cuando usábamos las distancias de Manhattan y con k-means usando las distancias euclídea y de Mahalanobis.

Using Mahalanobis Distance

Repetimos el truco de antes para transformar la matriz X (ya se mostró el código previamente, por lo que aquí se oculta; se puede ver completo en R Markdown, pero el objetivo es facilitar la lectura del report en HTML).

De nuevo utilizar la distancia de Mahalanobis parece una idea pésima, pues no se distinguen bien los grupos, no agrupa tan bien como el método de PAM usando las distancias Euclidea y de Manhattan. El mejor resultado de cluster para la técnica PAM, es por tanto cogiendo 3 clusters (matemáticamente serían 2, pero por nuestro objetivo del trabajo, como dijimos antes, nos encajan mejor 3), con la distancia de Manhattan.

Vamos a volver a visualizar estos resultados y a explicar cuáles son las caractrísticas de cada uno de estos clusters.

##    gdp_per_cap   inflation gni_per_cap     urban_p productivity credit_inf
## 33  -0.6435012 -0.70463209  -0.6525326 -0.10639363   -0.8740631 -1.2085410
## 83  -0.4445131  0.03126667  -0.4527898 -0.08911238   -0.5871169  0.7408338
## 56   1.3681889 -0.62749405   1.5045717  0.97800480    1.4412508  0.4159380
##      renew_en renew_electr electricity_2016 clean_fuels_2016 inf_mort_rate
## 33  1.5459515    1.3198344       -0.7815197       -1.1109502      1.668366
## 83 -0.5414891   -0.7114842        0.5907409        0.6797209     -0.442497
## 56 -0.6557915   -0.5388849        0.6553683        0.9318907     -0.924980
##    sanitation_services deaths_household_pollution      water  fert_rate
## 33          -1.1885475                  1.3631735 -1.1869961  1.4674467
## 83           0.3794512                 -0.6918295  0.3973235 -0.5448120
## 56           0.8566682                 -0.9854014  0.7934034 -0.6142002
##      life_exp urban_pop_growth    emp_agr    eco_agr vulner_empl
## 33 -1.8001742        0.8954380  1.5353881  0.8394029  1.34169151
## 83  0.5636106       -0.6947977 -0.3524337 -0.4231032  0.03917959
## 56  1.3556641       -0.7536954 -1.0060759 -0.8601245 -1.12229910
##       acc_own
## 33 -0.9453004
## 83 -0.8579625
## 56  1.2954307

Vemos que el primer cluster tiene un gdp per capita muy bajo (el más bajo de los 3 clusters), una inflación muy baja también, gni per capita muy bajo, población urbana baja pero cercana a la media, productividad muy baja (la más baja de todas), información al crédito muy muy baja, un altísimo uso de energías y de electricidad renovable, sin embargo bajísimo acceso a la electricidad y a combustibles limpios en el año 2016, una tasa de mortalidad infantil altísima, servicios sanitarios bajísimos, muy inferiores a la media, muertes por household pollution muy altas, bajíismo acceso al agua, altísima tasa de fertilidad, bajísima esperanza de vida en comparación con lo otros dos clusters; además, es el cluster que tiene una media de crecimiento de la población urbana más alta, con un altísimo peso del empleo y economía agrícola, un gran número de empleos vulnerables y un bajísimo porcentaje medio de account ownership. Vemos por lo tanto que este primer cluster representaría el tercer mundo, o aquellos países que peor posicionados están con respecto al efecto del spread. El cluster 2 vemos que tiene también unas variables económicas bajas, aunque no tan bajas como las del cluster 1… Muy similar al análisis que realizamos antes para el 3-mean custering con distancia euclídea; por lo tanto podríamos decir que este segundo cluster está formado por los países del segundo mundo, o aquellos países que no son por tanto ni los mejor posicionados con respecto al spread ni tampoco se están quedando atrás debido a este. EL último y tercer cluster está formado por los países que mejor puntuación media tienen en las variables económicas, así como en calidad de vida y sanidad (agua, servicios sanitarios…), como ya hemos dicho antes estos resultados son casi iguales que los que teníamos con el método de 3-mean clustering usando distancias euclídeas. Podemos identificar este tercer cluster con los países del primer mundo, o aquellos que estan mejor posicionados con respecto al spread.

Antes de pasar al siguiente método de clustering, vamos a dibujar en el mapa del mundo a los países según el grupo al que pertenecen por el mejor modelo de PAM clustering; es decir, 3 clusters utilizando la distancia de Manhattan.

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Vemos que por el método PAM los clusters salen muy parecidos, excepto alguna salvedad como Chile, que pasa de estar en el segundo mundo al primero.

Pruebas con Submuestras

Vamos a analizar cómo quedaría el mejor de los clusters por el método de PAM, es decir, por distancia Manhattan con 3 medioides, para las diferentes submuestras que se comentaron antes.

En primer lugar, vamos a ver cómo queda con los datos reducidos que sacamos en el estudio anterior de las variables:

## Silhouette of 182 units in 3 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        54        78        50 
## 0.1640951 0.3025651 0.3601817 
## Individual silhouette widths:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -0.5049  0.1665  0.3156  0.2773  0.4225  0.5761

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Quitando estas variables (ver arriba las variables concretas utilizadas en este primer subset, X_reduced) vemos que los resultados son curiosos; en primer lugar, Portugal y Grecia (especialmente llamativo el caso de esta última) estarían en el primer mundo, al igual que países como Arabia Saudí… El resto de países se comportan de forma parecida a cuando ejecutamos este método de clustering con los datos al completo.

Sólo con variables económicas:

## Silhouette of 182 units in 3 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        57        42        83 
## 0.2346880 0.4010943 0.3334542 
## Individual silhouette widths:
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
## -0.1777  0.2062  0.3485  0.3181  0.4418  0.6063

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Como ya comenté antes, el color y número de los clusters varía cada vez que se ejecuta el código, por lo que utilizamos referencias en países que tenemos claro a qué cluster deben pertenecer en función del análisis previo de los mismos (trabajo anterior cuyo link está expuesto en la introducción) para determinar las etiquetas que les ponemos a estos clusters.

Según el método de clustering utilizando estas variables, los países que forman el primer mundo son prácticamente los mismos que para el subset de datos anterior; Grecia, por ejemplo, sí que se mueve de cluster y entraría en el segundo mundo. La mayor parte de África y algunos países de oriente y Asia como Afghanistán pertenecen al tercer mundo según las variables puramente económicas. Probemos ahora con las variables relacionadas sólo con la salud:

## Silhouette of 182 units in 3 clusters from silhouette.default(x = groups, dist = d) :
##  Cluster sizes and average silhouette widths:
##        55        52        75 
## 0.1699682 0.2425302 0.1790369 
## Individual silhouette widths:
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
## -0.45036  0.07378  0.26223  0.19444  0.35968  0.51557

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Según las variables de salud, vemos unos resultados bastante similares en algunos aspectos a los que observábamos con el método de K-Means con 3 centroides: las diferencias entre los grupos se disipan un poco, en el sentido de que las medias para éstas variables de un cluster a otro no difieren tanto. En especial debe notarse una mayor inclusión de países en el cluster de los países del primer mundo, en el que entran ahora Argentina, Chile, Turquía, y otros países que no estaban incluidos en este cluster utilizando todas las variables, la matriz reducida o la matriz de variables económicas. Esto es una buena noticia, pues significa que, pese a que siga habiendo millones de personas viviendo con pocos nutrientes, alta contaminación que provoca muertes, con un mal acceso al agua y a la sanidad etc, cada vez más países gozan de estos derechos humanos.

5. Model-Based Clustering (Gaussian Mixtures)

Primero vamos a dejar que sea el propio modelo el que optimice el número de clusters (anticipamos que seguramente lo divida en 2, como hemos visto antes parece que al dividirlos en 2 matemáticamente es con el número de clusters que más tienen en común los países de cada cluster, teóricamente); luego después ajustaremos un modelo de 3 clusters, por los motivos que ya explicamos arriba, pues es con el número de grupos que más cómodos nos sentimos.

require(mclust)
fit.m <- Mclust(X)
summary(fit.m) # seems 2 clusters fit well
## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust VVE (ellipsoidal, equal orientation) model with 2 components: 
## 
##  log.likelihood   n  df       BIC       ICL
##       -1952.266 182 295 -5439.714 -5440.387
## 
## Clustering table:
##   1   2 
## 106  76
# we can also specify tghe number of groups
# fit.m <- Mclust(X, G=3)

plot(fit.m, what = "classification", main="")

plot(fit.m, what = "BIC", main="")

# By default, the models considered are:
# "EII": spherical, equal volume 
# "VII": spherical, unequal volume 
# "EEI": diagonal, equal volume and shape
# "VEI": diagonal, varying volume, equal shape
# "EVI": diagonal, equal volume, varying shape 
# "VVI": diagonal, varying volume and shape 
# "EEE": ellipsoidal, equal volume, shape, and orientation 
# "EEV": ellipsoidal, equal volume and equal shape
# "VEV": ellipsoidal, equal shape 
# "VVV": ellipsoidal, varying volume, shape, and orientation

groups = fit.m$classification

centers = fit.m$parameters$mean

i=1  # plottinng the centers in cluster 1
barplot(centers[,i], las=2, col="darkblue", main=paste("Cluster", i))

# insights?

i=2  # plottinng the centers in cluster 2
barplot(centers[,i], las=2, col="darkblue", main=paste("Cluster", i))

# scatter plot:
fviz_cluster(fit.m,geom = c("point"), ellipse.type = 'norm',pointsize=1) +
  theme_minimal()+geom_text(label=countries,hjust=0, vjust=0,size=2,check_overlap = T)+scale_fill_brewer(palette="Paired")

Vemos que la agrupación que hace este método no se diferencia mucho a la que probamos con métodos anteriores cogiendo 2 clusters.

Ahora vamos a ajustar un modelo de 3 clusters:

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust VEV (ellipsoidal, equal shape) model with 3 components: 
## 
##  log.likelihood   n  df      BIC       ICL
##       -1184.192 182 718 -6104.86 -6104.961
## 
## Clustering table:
##  1  2  3 
## 73 39 70

Analizaremos ahora las características de los clusters basándonos en los 3 gráficos de barras que hemos dibujado arriba. En ellos se ve la media (recordemos que todas las variables están normalizadas) que tiene cada uno de los clusters en las 21 variables sobre los países que estamos considerando.

El primer cluster vemos que tiene unas variables económicas (gdp per capita, gni per capital, productividad, información crediticia…) bastante bajas; por debajo de la media. Llaman la atención el bajísimo acceso a la electricidad y a combustibles limpios en el año 2016, la alta tasa de mortalidad infantil, el bajo acceso a servicios sanitarios y al agua, la alta tasa de fertilidad… En general, vemos de nuevo que este primer cluster correspondería a aquellos países que pertenecen al tercer mundo, es decir lo forman los países peor posicionados con respecto al spread.

El segundo cluster está claro que pertenece a los países que forman parte del primer mundo, es decir aquellos países que están mejor posicionados con respecto al spread. Esto lo vemos claro en el gráfico de barras, pues es el grupo que mayor media tiene en las variables económicas y otras variables que miden la calidad de vida, mientras que tiene una tasa bajísima de mortalidad infantil, muy pocos empleos vulnerables etc.

Por último, el tercer y último cluster está entre medias de los otros dos en casi todas las variables, podemos ver que no puntúa tan alto en las variables “buenas” (aquellas asociadas con buen posicionamiento con respecto al spread) como el segundo cluster, ni tan mal como el primero; por otro lado, no puntúa tan bajo en las variables “malas” como el segundo, ni tan alto como el primero. Por este motivo esto nos lleva a decir que este tercer cluster son los países que están posicionados de forma intermedia con respecto al spread.

Ahora visualizamos los resultados de este último modelo de cluster en el mapa del mundo, para ver si existen diferencias significativas con respecto al mejor modelo de PAM y de k-means. Primero cogeremos los grupos, y cambiaremos el orden de estos para que visualmente sea más fácilmente comparable con los otros dos gráficos que dibujamos antes, pues recordemos que en este caso el primer mundo por ejemplo en lugar de ser el cluster 3 es el 2.

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Vemos que en este caso los resultados también se asemejan a los obtenidos por los otros dos métodos de clustering; se asemeja más a la separación de grupos que había hecho el k-means. Sin embargo, este modelo es algo más fino y países como India los mete en el cluster 1, es decir en los países del tercer mundo.

Pruebas con Submuestras

Primero probaremos con la submuestra reducida:

fit.m <- Mclust(X_reduced, G=3)
summary(fit.m) 
## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust VEV (ellipsoidal, equal shape) model with 3 components: 
## 
##  log.likelihood   n  df       BIC       ICL
##       -406.7847 182 428 -3040.884 -3041.013
## 
## Clustering table:
##  1  2  3 
## 74 69 39
plot(fit.m, what = "classification", main="")

plot(fit.m, what = "BIC", main="")

groups = fit.m$classification

centers = fit.m$parameters$mean

i=1  # plottinng the centers in cluster 1
barplot(centers[,i], las=2, col="darkblue", main=paste("Cluster", i))

# insights?

i=2  # plottinng the centers in cluster 2
barplot(centers[,i], las=2, col="darkblue", main=paste("Cluster", i))

# insights?

i = 3
barplot(centers[,i], las=2, col="darkblue", main=paste("Cluster", i))

# scatter plot:
fviz_cluster(fit.m,geom = c("point"), ellipse.type = 'norm',pointsize=1) +
  theme_minimal()+geom_text(label=countries,hjust=0, vjust=0,size=2,check_overlap = T)+scale_fill_brewer(palette="Paired")

matched <- joinCountryData2Map(data.frame(country=countries, value=groups), joinCode="NAME", nameJoinColumn="country")
## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data
mapCountryData(matched, nameColumnToPlot="value", mapTitle="Clusters for Countries in the World using Model-Based Clustering", 
               catMethod = "categorical", colourPalette = "topo")

En este caso no hay demasiado que comentar, salvo, tal vez, que según este método de clustering, utilizando los datos reducidos, países como Colombia caerían al cluster de los países del tercer mundo. Este es también el caso de Corea del Sur… es curioso, pues por los métodos anteriores con estos datos no se observó un comportamiento similar.

Sólo con los datos económicos:

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust EVE (ellipsoidal, equal volume and orientation) model with 3
## components: 
## 
##  log.likelihood   n df       BIC       ICL
##       -643.7723 182 90 -1755.905 -1766.705
## 
## Clustering table:
##  1  2  3 
## 90 58 34

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Arabia Saudí, que antes formaba parte del primer mundo con datos económicos, ahora formaría parte del segundo mundo, al igual que la mayor parte del golfo pérsico. Algunos países como Colombia o Ecuador caerían en el primer mundo en términos puramente económicos.

Sólo con datos de salud:

## ---------------------------------------------------- 
## Gaussian finite mixture model fitted by EM algorithm 
## ---------------------------------------------------- 
## 
## Mclust VEV (ellipsoidal, equal shape) model with 3 components: 
## 
##  log.likelihood   n  df       BIC      ICL
##       -670.8148 182 179 -2273.147 -2274.68
## 
## Clustering table:
##  1  2  3 
## 59 42 81

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

COn los datos únicamente de salud podemos ver que la mayoría entran dentro del primer mundo si sólo atendemos a estas variables. Nos llama especialmente la atención que algunos de los países que mejor puntuaban tanto en el factor spread general como en el factor spread económico, como Noruega, Suecia, Islandia… se posicionarían en el segundo mundo según el Model-Based Clustering con los datos de salud. Esto es raro, pues quedarían en el mismo cluster que países como India…

Por último, antes de comparar todos los métodos de cluster y de quedarnos con el modelo (o los modelos) que más nos convenza, vamos a implementar Hierarchical Clustering.

6. Hierarchical Clustering.

En el caso del hierarchical clustering, y dado que ya hemos comentado sobradamente que 3 clusters nos encajan mejor que 2 para el objetivo de este estudio, vamos a probar solamente haciendo 3 clusters (cortando el árbol de tal forma que dejemos 3 grupos, ya que en este método de cluster, en función de donde “cortes el árbol” obtienes un número u otro de clusters).

Hierarchical Clustering with Euclidean Distance

if(!require(ape)) {
  
  install.packages("ape")
  
  require(ape)
  
}
## Loading required package: ape
## Warning: package 'ape' was built under R version 3.5.2
d <- dist(X, method = "euclidean") # distance matrix using euclidean distance
# we can choose any of "euclidean", "maximum", "manhattan", "canberra", "binary" or "minkowski" 


fit.h <- hclust(d, method="ward.D2") # cluster using Ward linkage (based on minimizing variance)

dendro <- as.dendrogram(fit.h)

fviz_dend(fit.h, k=3, rect = TRUE) # dendrogam

groups.h <- cutree(fit.h, k=3) 

X_grouped <- X

X_grouped$group <- groups.h

En este caso dibujamos el mapa del mundo con los clusters para cada modelo de Hierarchical Clustering, pues es la forma más clara de determinar, visualmente, si nos encajan los resultados obtenidos en cada caso.

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Vemos que los resultados son muy similares a los obtenidos por el método de k-means y PAM anteriormente (más parecido al primero).

Haremos la media de las variables para cada uno de los clusters para poder determinar con algo más de claridad las características de cada uno de ellos.

## # A tibble: 3 x 22
##   group gdp_per_cap inflation gni_per_cap urban_p productivity credit_inf
##   <int>       <dbl>     <dbl>       <dbl>   <dbl>        <dbl>      <dbl>
## 1     1      -0.634    0.405       -0.652  -0.873      -0.832      -0.758
## 2     2      -0.253   -0.0575      -0.255   0.110      -0.0769      0.266
## 3     3       1.86    -0.555        1.90    1.24        1.71        0.578
## # ... with 15 more variables: renew_en <dbl>, renew_electr <dbl>,
## #   electricity_2016 <dbl>, clean_fuels_2016 <dbl>, inf_mort_rate <dbl>,
## #   sanitation_services <dbl>, deaths_household_pollution <dbl>,
## #   water <dbl>, fert_rate <dbl>, life_exp <dbl>, urban_pop_growth <dbl>,
## #   emp_agr <dbl>, eco_agr <dbl>, vulner_empl <dbl>, acc_own <dbl>

Vemos que el primer grupo son los países del tercer mundo, con un gdp per capita y un gni per capita muy bajos, con poquísima población urbana (siempre hablamos en medias y con respecto a los otros grupos), con bastante inflación (mucha más que los otros dos grupos); además es el grupo con menor información de crédito, y menor acceso a electricidad y a combustibles limpios en 2016. Curiosamente los países del tercer mundo son los que mayor uso de electricidad y energías renovables tienen. Sin embargo, tienen una mortalidad infantil altísima, unos servicios sanitarios muy bajos, así como bajísimo acceso al agua, muchas muertes por household pollution, una esperanza de vida muy inferior a la de los otros grupos, un gran crecimiento de la población urbana… El empleo y la economía agrícola tienen un gran peso en las sociedades que forman parte de este primer cluster, que tienen además un gran número de empleos vulnerables y un bajo porcentaje de adultos con cuenta bancaria en su posesión.

El segundo grupo serían los países que forman el segundo mundo, ya que vemos que tienen las características “buenas” (sanidad, agua, variables económicas, posesión de cuenta bancaria, esperanza de vida) más altas que el cluster 1 pero más bajas que el cluster 2. Además, las características “negativas” (muertes por contaminación, empleo y economía agrícolas, empleo vulnerable, inflación…) son más bajas en este segundo cluster que en el primero, pero más altas que en el tercero.

Por último, el tercer y último cluster lo forman aquellos países que tienen mayor puntuación en las variables de índole positiva (productividad, esperanza de vida… todas las mencionadas antes), mientras que tienen la puntuación más baja en las variables de índole negativa, como mortalidad infantil o empleo vulnerable (además de las ya mencionadas previamente como negativas). Podemos concluir por tanto que este tercer cluster lo forman los países del primer mundo, es decir aquellos que mejor posicionados están con respecto al spread.

Anticipamos que este es el modelo de Hierarchical Clustering que más nos encaja tras realizar todo el resto del análisis que se puede encontrar más abajo. Por este motivo probaremos este modelo de clustering sobre los 3 subconjuntos de datos.

Pruebas con Subconjuntos

Con los datos reducidos:

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Con los datos reducidos, en este caso, hay varias cosas que nos llaman la atención. En primer lugar, ahora ya no es sólo Portugal, sino que España e Italia se suman a los países que pasarían a ser parte del segundo mundo según este subconjunto de variables. Vemos, por lo tanto, que algunas de las variables eliminadas (verlas arriba) eran clave para que estos países formen parte del primer mundo según el Hierarchical Clustering con 3 clusters.

Probamos ahora con datos económicos:

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

La mayoría de lo que dijimos para los datos reducidos es cierto también para el Hierarchical CLustering con 3 clusters usando sólo los datos económicos. Vemos que países como España, Portugal e Italia pasarían a formar parte del segundo mundo, quedando por detrás de países como Arabia Saudí, que pasaría a formar parte del primer mundo.

Por último, realizamos estas pruebas sobre el subconjunto de datos de salud:

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Es bastante similar a los que hemos visto por otros métodos para este subconjunto de datos. Podemos ver que la mayoría de los países formarían parte del primer mundo, pero hay algunos aspectos de este hierarchical clustering que no nos cuadran mucho… Canadá, Noruega, Islandia… forman parte del segundo mundo en términos de salud utilizando este método de clustering. Esto no es lo esperado, ya que son países que, además, tenían de las mejores puntuaciones en los factores (las variables latentes existentes en los datos completos), que representaban el posicionamiento con respecto al spread; sin embargo, parece que cuando aislamos las variables económicas de los datos y dejamos solamente las de salud, no tienen tan buena puntuación como otros países como Rusia, España, Argentina o Chile. Por puntuación, entiéndase que no tienen medias tan altas para las variables “buenas” (acceso al agua, acceso a la sanidad, esperanza de vida…) ni medias tan poco bajas para las variables “malas” (mortalidad infantil, muerte por contaminación de hogares…).

Ahora probaremos cogiendo otras distancias, en este caso la distancia de Manhattan.

Hierarchical Clustering Using Manhattan Distance

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Utilizando las distancias de Manhattan nos convence menos el hierarchical cluster, pues como vemos algunos países como Rusia o Argentina (el uno con un régimen político de dudosa legitimidad moral y el otro con problemas económicos cíclicos) caerían dentro de los países del primer mundo, lo cual no nos cuadra. Probaremos ahora con otro tipo de distancia.

Hierarchical Clustering Using Maximum Distance.

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Vemos que en este caso, utilizando como medida de distancias “Maximum”, separa los países parecido a los otros métodos utilizados, con algunas salvedades como que más países africanos caerían dentro del segundo mundo, y más países del primer mundo también quedarían dentro de este segundo mundo.

Hierarchical Clustering Using Minkowski Distance

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Vemos que apenas cambia nada (salvo algunos países como India) con respecto al Hierarchical Clustering utilizando las distancias máximas.

6 (bis). Advanced Hierarchical Clustering.

Using Euclidean Distance.

require(factoextra)

fit.hc <- eclust(X, 'hclust', k=3, hc_metric = 'euclidean', hc_method = "ward.D2", graph=F)
groups = fit.hc$cluster
fviz_dend(fit.hc, rect = TRUE) # dendrogam

fviz_silhouette(fit.hc)
##   cluster size ave.sil.width
## 1       1   57          0.30
## 2       2   93          0.27
## 3       3   32          0.44

matched <- joinCountryData2Map(data.frame(country=countries, value=groups), joinCode="NAME", nameJoinColumn="country")
## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data
mapCountryData(matched, nameColumnToPlot="value", mapTitle="Clusters for Countries in the World using Hierarchical Clustering", 
               catMethod = "categorical", colourPalette = "topo")

Los resultados de este cluster nos salen muy parecidos a los que obtuvimos antes, como podemos ver en el mapa del mundo.

Using Manhattan Distance.

##   cluster size ave.sil.width
## 1       1   56          0.34
## 2       2   65          0.26
## 3       3   61          0.37

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

El average silhouette width sale algo mejor que en el anterior caso, aunque solo 0.01 más; sin embargo vemos en el mapa que nos cuadra menos con la idea que tenemos, pues pasa algo parecido a lo que vimos antes, y es que países como Rusia y Argentina entran dentro del cluster 3 (primer mundo), lo cual no es muy consistente con las ideas que tenemos de estos países.

Using Canberra Distance.

##   cluster size ave.sil.width
## 1       1   49          0.47
## 2       2   49         -0.03
## 3       3   84          0.41

## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data

Utilizando la distancia de Canberra los resultados son aún menos coherentes, pues muchos países que según este modelo de clustering entrarían en el primer mundo (como Brasil), está claro que no pueden formar parte del mismo grupo que países como Islandia, Singapur, Luxemburgo…

Using the “Complete” Method instead of Ward.D2.

Ahora vamos a probar, ya como último modelo de clustering (antes de comparar matemáticamente, con algunas estadísticas, todos los modelos de cluster que nos han convencido), el hierarchical clustering utilizando distancias euclídeas pero con el método de Complete.

fit.hc <- eclust(X, 'hclust', k=3, hc_metric = 'euclidean', hc_method = "complete", graph=F)
groups = fit.hc$cluster
fviz_dend(fit.hc, rect = TRUE) # dendrogam

fviz_silhouette(fit.hc)
##   cluster size ave.sil.width
## 1       1   65          0.28
## 2       2  109          0.28
## 3       3    8          0.43

matched <- joinCountryData2Map(data.frame(country=countries, value=groups), joinCode="NAME", nameJoinColumn="country")
## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data
mapCountryData(matched, nameColumnToPlot="value", mapTitle="Clusters for Countries in the World using Hierarchical Clustering", 
               catMethod = "categorical", colourPalette = "topo")

Vemos que con este método se quedan solos en el cluster 3 (primer mundo) los países con variables económicas muy por encima del resto, teniendo menos peso a la hora de formar parte de ese cluster el resto de variables del nivel de vida que comentamos antes.

7. Comparing The Best Clustering Models.

Como para sacar estas estadísticas necesitamos utilizar una única matriz de distancias para todos los clusters, y la matriz de distancias que mejor resultados nos ha dado en general ha sido utilizando la distancia euclídea (excepto en el caso de PAM que nos cuadraban mejor los resultados con la distancia de Manhattan), será esta la distancia utilizada para sacar la matriz de distancias.

Con la función cluster.stats comparamos el acuerdo que hay entre 2 composiciones de cluster. Con el corrected Rand Index medimos este nivel de acuerdo, de -1 (no acuerdo en absoluto) a 1 (totalmente de acuerdo). Por otro lado, cuanto más cerca esté el Meila´s VI del 0, más similares son los clusters.

fit.k_means = kmeans(X, centers=3, nstart=100)

groups_kmeans <- fit.k_means$cluster

fit.pam <- pam(X, k=3, metric = "euclidean")

groups_pam <- fit.pam$clustering

fit.m_clust <- Mclust(X, G=3)

groups_mclust <- fit.m$classification

d <- dist(X, method = "euclidean")

fit.h <- hclust(d, method="ward.D2")

groups.h <- cutree(fit.h, k=3)

Comparación K-Means vs PAM

if(!require(fpc)) {
  
  install.packages("fpc")
  
  library(fpc)
}


clust_stats = cluster.stats(d, groups_kmeans, groups_pam)

# Among the values returned by previous function, there are two indexes to assess the similarity of two clustering:
#      Corrected Rand index and VI Score
# Rand Index should be maximized and VI score should be minimized
clust_stats$corrected.rand
## [1] 0.5737718
clust_stats$vi
## [1] 0.7556367

Para poder poner en perspectiva los resultados mostrados arriba, es importante comparar todo el resto de clusters con todos.

Debajo iremos sacando, para los distintos clusters, el corrected rand index y el VI score de arriba. Como el código es igual al que viene aquí arriba, no lo pondremos y se enseñará sólo el output.

Comparación K-Means vs Hierarchical Clustering

## [1] 0.8568725
## [1] 0.3373642

Vemos que el K-Means y el Hierarchical Clustering se parecen mucho, pues tienen un Corrected Rand Index de 0.86, bastante por encima del 0.57 del K-means y PAM; el Meila´s VI es bastante inferior a este anterior caso. Podemos decir, por tanto, que el K-Means y el Hierarchical Clustering tienen bastante similaridad, agrupan a la mayoría de los países en el mismo cluster.

Comparación PAM vs Hierarchical Clustering

## [1] 0.4903959
## [1] 0.8554803

Sin embargo, PAM y Hierarchical Clustering no se parecen demasiado, menos que K-Means y PAM y mucho menos que K-Means y Hierarchical Clustering. De momento parece que entre K-Means y PAM nos quedaríamos con el primero; los resultados de ambos debieran ser muy parecidos, con la diferencia de que PAM encuentra sus medioides entre los países que existen, mientras que los centroides de K-Means son ficticios (no son observaciones reales); esto puede significar una mejor interpretabilidad de los clusters de PAM, sin embargo puede llevar (especialmente al tratar con datos de países, pues hay pocas observaciones y los medioides no tienen tanta libertad para elegirse) a una peor agrupación de los países.

Comparación Model-Based Clustering vs. Hierarchical Clustering.

## [1] 0.3828165
## [1] 1.139037

Vemos que estos dos tipos de cluster se parecen algo más que PAM y Hierarchical Clustering y que PAM y K-Means, sin embargo se parecen bastante menos que Hierarchical Clustering y K-Means.

COmparación K-Means vs. Model-Based Clustering.

## [1] 0.3963289
## [1] 1.093381

Vemos que la similaridad de estos dos clusters es muy similar a la de K-Means y PAM; es bastante inferior a la de K-Means y Hierarchical Clustering.

Comparación PAM vs. Model-Based Clustering

## [1] 0.4160893
## [1] 1.121986

Como vemos, PAM no tiene mucho que ver con Model-Based Clustering, en términos de los resultados de la agrupación, pues estas dos agrupaciones son bastante diferentes, a juzgar por el corrected rand index y el Meila´s VI score. Podemos determinar que PAM no tiene demasiado en común con ninguno de los otros métodos de Clustering.

8. Conclusiones

Finalmente, se concluye que los dos modelos de Clustering que más tienen en común son K-Means y Hierarchical Clustering; serán por lo tanto los dos modelos de clustering con los que nos quedamos. Volveremos a graficar sus respectivos resultados sobre el mapa del mundo para visualizar las diferencias:

matched1 <- joinCountryData2Map(data.frame(country=countries, value=groups.h), joinCode="NAME", nameJoinColumn="country")
## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data
mapCountryData(matched1, nameColumnToPlot="value", mapTitle="Final Clusters for Countries in the World using Hierarchical Clustering", 
               catMethod = "categorical", colourPalette = "topo")

matched2 <- joinCountryData2Map(data.frame(country = countries, value = groups_kmeans), joinCode = "NAME",
                                nameJoinColumn = "country")
## 176 codes from your data successfully matched countries in the map
## 6 codes from your data failed to match with a country code in the map
## 67 codes from the map weren't represented in your data
mapCountryData(matched2, nameColumnToPlot="value", mapTitle="Final Clusters for Countries in the World using K-Means", 
               catMethod = "categorical", colourPalette = "topo")

Vemos que las diferencias entre ambos métodos de cluster son poco significativas. Algunas que nos llaman la atención por ejemplo son:

Vemos, por tanto, que pese a ser similares, el Hierarchical Clustering es algo más exigente a la hora de incluir un país en un cluster avanzado, y tiende a meter a los países en clusters más bajos que el K-Means. En las discrepancias que vemos, nos decidimos quedar con el Hierarchical Clustering, pues nos sentimos más cómodos incluyendo a países como Grecia o Eslovenia en el cluster 2 que en el cluster 1 (no son los países que mejor posicionados están con respecto al spread, y tienen gran diferencia, especialmente en los últimos años – piense el lector en el caso de Grecia) con respecto a los países más desarrollados como Islandia, Singapur o Luxemburgo. Esto es coherente pues generalmente el Hierarchical Clustering es algo más fino que el K-Means clustering; es computacionalmente más costoso, pero esto con los datos que tratamos en este caso no es un problema, pues el número de países existentes es relativamente limitado. En algunos casos nos surgen dudas, como es el caso de Estonia, pues es un país que en los últimos años ha llevado a cabo reformas de transformación digital; además es el país con mayor número de puestos de carga de coches eléctricos de Europa. Sin embargo, históricamente no ha tenido una calidad de vida como la que han gozado los habitantes de países como Luxemburgo. En cualquier caso, los indicadores de calidad de vida como fertility rate o esperanza de vida, así como el acceso a una sanidad de calidad, tardan varios años en estabilizarse, y transformaciones radicales como la observada en Singapur son sucesos remotos. Por este motivo, consideramos que de momento Estonia pertenece a los países del segundo mundo, aquellos que están bien posicionados con respecto al spread en comparación con los países menos desarrollados (cluster 1 – ver gráficos arriba), pero a los que aún les falta para llegar a posicionarse como los países del primer mundo (cluster 3).

Las descripciones de las características de cada cluster (de los dos de arriba) están descritas en sus correspondientes apartados.

Antes de finalizar, sacaremos algunas estadísticas para validar los clusters, con la ayuda de la librería clValid.

if(!require(clValid)) {
  
  install.packages("clValid")
  
  library(clValid)
}

Primero veremos la estabilidad del hierarchical cluster cuando eliminamos cada una de las columnas. Para ello, vamos a definir primero un poco las medidas de estabilidades que se incluyen en este paquete, así como su significado, para hacer los resultados más fáciles de interpretar.

Extraigo un texto de este link en el que se explica un poco en qué consisten estas medidas.

The cluster stability measures includes:

The average proportion of non-overlap (APN) The average distance (AD) The average distance between means (ADM) The figure of merit (FOM)

The APN, AD, and ADM are all based on the cross-classification table of the original clustering with the clustering based on the removal of one column.

The APN measures the average proportion of observations not placed in the same cluster by clustering based on the full data and clustering based on the data with a single column removed.

The AD measures the average distance between observations placed in the same cluster under both cases (full dataset and removal of one column).

The ADM measures the average distance between cluster centers for observations placed in the same cluster under both cases.

The FOM measures the average intra-cluster variance of the deleted column, where the clustering is based on the remaining (undeleted) columns. It also has a value between zero and 1, and again smaller values are preferred.

The values of APN, ADM and FOM ranges from 0 to 1, with smaller value corresponding with highly consistent clustering results. AD has a value between 0 and infinity, and smaller values are also preferred.

stab <- matrix(0,nrow=ncol(X),ncol=4)
colnames(stab) <- c("APN","AD","ADM","FOM")

clusterObj <- hclust(d, method="ward.D2")

cluster <- cutree(clusterObj,3)


for (del in 1:ncol(X)) {
  matDel <- X[,-del]               
  DistDel <- dist(matDel,method="euclidean")
  clusterObjDel <- hclust(DistDel, method="ward.D2")
  clusterDel <- cutree(clusterObjDel,3)
  stab[del,] <- stability(X, DistDel, del, cluster, clusterDel)
}
##              APN       AD        ADM       FOM
##  [1,] 0.00000000 3.845971 0.00000000 0.4790769
##  [2,] 0.18072476 3.931146 0.84086964 0.9457717
##  [3,] 0.00000000 3.851354 0.00000000 0.4421126
##  [4,] 0.05138175 3.859479 0.28983190 0.7169767
##  [5,] 0.03147525 3.874390 0.12078494 0.5508391
##  [6,] 0.07826347 3.876759 0.43917373 0.8640439
##  [7,] 0.02120686 3.821570 0.15418012 0.7298832
##  [8,] 0.20311946 3.958159 0.96000956 0.9800631
##  [9,] 0.00000000 3.845271 0.00000000 0.5449667
## [10,] 0.16613494 4.045979 0.74384031 0.4679173
## [11,] 0.04227147 3.917577 0.18755639 0.5726548
## [12,] 0.20223014 4.109744 0.91756607 0.4669631
## [13,] 0.05723074 3.919855 0.31150259 0.6000425
## [14,] 0.01079622 3.865224 0.08060839 0.5747106
## [15,] 0.00000000 3.823502 0.00000000 0.6095749
## [16,] 0.16613494 4.033320 0.74464163 0.5087956
## [17,] 0.09913742 3.960724 0.61607958 0.8022540
## [18,] 0.00000000 3.834813 0.00000000 0.5780768
## [19,] 0.01064560 3.812588 0.05448244 0.7184186
## [20,] 0.02144183 3.863302 0.11143867 0.5771236
## [21,] 0.17251566 3.981086 0.77102291 0.6454997
print(paste("el mínimo AD es quitando la columna" , which(stab[,1] == min(stab[,1])), "con un valor de ", stab[ which(stab[,1] == min(stab[,1])), 1]))
## [1] "el mínimo AD es quitando la columna 1 con un valor de  0" 
## [2] "el mínimo AD es quitando la columna 3 con un valor de  0" 
## [3] "el mínimo AD es quitando la columna 9 con un valor de  0" 
## [4] "el mínimo AD es quitando la columna 15 con un valor de  0"
## [5] "el mínimo AD es quitando la columna 18 con un valor de  0"
print(paste("el mínimo APN es quitando la columna" , which(stab[,2] == min(stab[,2])), "con un valor de ", stab[which(stab[, 2] == min(stab[,2])), 2]))
## [1] "el mínimo APN es quitando la columna 19 con un valor de  3.81258817255093"
print(paste("el mÁXIMO ADM es quitando la columna" , which(stab[, 3] == max(stab[, 3])), "con un valor de ", stab[which(stab[, 3] == max(stab[, 3])), 3]))
## [1] "el mÁXIMO ADM es quitando la columna 8 con un valor de  0.960009555914833"

EN este output vemos cuáles son los valores mínimos de las primeras dos métricas y el mayor de la tercera. Los valores más bajos en los dos primeros casos corresponden a las columnas que menos afectan a la estabilidad del cluster al estar ausentes. En el tercer caso, el valor máximo que aparece corresponde a la columna que más afecta (más desestabiliza la estructura de los clusters).

Ahora realizaremos un análisis parecido con el k-means.

stab_k <- matrix(0,nrow=ncol(X),ncol=4)
colnames(stab_k) <- c("APN","AD","ADM","FOM")

clusterObj <- kmeans(X, 3, iter.max = 100, nstart = 100)

cluster = clusterObj$cluster

for (del in 1:ncol(X)) {
  matDel <- X[,-del]               
  clusterObjDel <- kmeans(X, 3, iter.max = 100, nstart = 100)
  DistDel <- dist(matDel, "euclidean")
  clusterDel <- clusterObjDel$cluster
  stab_k[del,] <- stability(X, DistDel, del, cluster, clusterDel)
}
##       APN       AD ADM       FOM
##  [1,]   0 3.834387   0 0.5134857
##  [2,]   0 3.706360   0 0.9438560
##  [3,]   0 3.838076   0 0.4814029
##  [4,]   0 3.769938   0 0.7141293
##  [5,]   0 3.822644   0 0.5508920
##  [6,]   0 3.723421   0 0.8465519
##  [7,]   0 3.779840   0 0.7294605
##  [8,]   0 3.662082   0 0.9792500
##  [9,]   0 3.845840   0 0.5164883
## [10,]   0 3.848425   0 0.4785885
## [11,]   0 3.843286   0 0.5198168
## [12,]   0 3.858936   0 0.4375166
## [13,]   0 3.813313   0 0.6079517
## [14,]   0 3.842742   0 0.5291888
## [15,]   0 3.826271   0 0.5769619
## [16,]   0 3.845963   0 0.4878497
## [17,]   0 3.757791   0 0.7566231
## [18,]   0 3.828760   0 0.5875240
## [19,]   0 3.796190   0 0.7106297
## [20,]   0 3.825697   0 0.5745140
## [21,]   0 3.798726   0 0.6429687
##  [1] "el mínimo AD es quitando la columna 1 con un valor de  0" 
##  [2] "el mínimo AD es quitando la columna 2 con un valor de  0" 
##  [3] "el mínimo AD es quitando la columna 3 con un valor de  0" 
##  [4] "el mínimo AD es quitando la columna 4 con un valor de  0" 
##  [5] "el mínimo AD es quitando la columna 5 con un valor de  0" 
##  [6] "el mínimo AD es quitando la columna 6 con un valor de  0" 
##  [7] "el mínimo AD es quitando la columna 7 con un valor de  0" 
##  [8] "el mínimo AD es quitando la columna 8 con un valor de  0" 
##  [9] "el mínimo AD es quitando la columna 9 con un valor de  0" 
## [10] "el mínimo AD es quitando la columna 10 con un valor de  0"
## [11] "el mínimo AD es quitando la columna 11 con un valor de  0"
## [12] "el mínimo AD es quitando la columna 12 con un valor de  0"
## [13] "el mínimo AD es quitando la columna 13 con un valor de  0"
## [14] "el mínimo AD es quitando la columna 14 con un valor de  0"
## [15] "el mínimo AD es quitando la columna 15 con un valor de  0"
## [16] "el mínimo AD es quitando la columna 16 con un valor de  0"
## [17] "el mínimo AD es quitando la columna 17 con un valor de  0"
## [18] "el mínimo AD es quitando la columna 18 con un valor de  0"
## [19] "el mínimo AD es quitando la columna 19 con un valor de  0"
## [20] "el mínimo AD es quitando la columna 20 con un valor de  0"
## [21] "el mínimo AD es quitando la columna 21 con un valor de  0"
## [1] "el mínimo APN es quitando la columna 8 con un valor de  3.66208195289207"
##  [1] "el mÁXIMO ADM es quitando la columna 1 con un valor de  0" 
##  [2] "el mÁXIMO ADM es quitando la columna 2 con un valor de  0" 
##  [3] "el mÁXIMO ADM es quitando la columna 3 con un valor de  0" 
##  [4] "el mÁXIMO ADM es quitando la columna 4 con un valor de  0" 
##  [5] "el mÁXIMO ADM es quitando la columna 5 con un valor de  0" 
##  [6] "el mÁXIMO ADM es quitando la columna 6 con un valor de  0" 
##  [7] "el mÁXIMO ADM es quitando la columna 7 con un valor de  0" 
##  [8] "el mÁXIMO ADM es quitando la columna 8 con un valor de  0" 
##  [9] "el mÁXIMO ADM es quitando la columna 9 con un valor de  0" 
## [10] "el mÁXIMO ADM es quitando la columna 10 con un valor de  0"
## [11] "el mÁXIMO ADM es quitando la columna 11 con un valor de  0"
## [12] "el mÁXIMO ADM es quitando la columna 12 con un valor de  0"
## [13] "el mÁXIMO ADM es quitando la columna 13 con un valor de  0"
## [14] "el mÁXIMO ADM es quitando la columna 14 con un valor de  0"
## [15] "el mÁXIMO ADM es quitando la columna 15 con un valor de  0"
## [16] "el mÁXIMO ADM es quitando la columna 16 con un valor de  0"
## [17] "el mÁXIMO ADM es quitando la columna 17 con un valor de  0"
## [18] "el mÁXIMO ADM es quitando la columna 18 con un valor de  0"
## [19] "el mÁXIMO ADM es quitando la columna 19 con un valor de  0"
## [20] "el mÁXIMO ADM es quitando la columna 20 con un valor de  0"
## [21] "el mÁXIMO ADM es quitando la columna 21 con un valor de  0"

En este caso vemos un suceso muy curioso, y es que el AD y el ADM es de 0 quitando muchas columnas, lo cual nos viene a decir que los 3 clusters por k-means son muy estables, pues quitando cualquiera de estas columnas que vienen arriba, los clusters seguirían igual de estables. Además, vemos que el ADM tiene un máximo de 0 en todas las columnas. COnsiderando que esta medida puede ir de 0 a infinito, es todo un logro para el cluster de k-means tener una estabilidad perfecta según esta medida, dando igual las variables que se eliminen del modelo. El mínimo APN es menor que para el hierarchical clustering con 3 clusters.

Podemos determinar, por lo tanto, que el cluster más estable es el k-means cluster con 3 centroides.